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Abstract 

Os '• 

■ The shape dynamics of fluid vesicles is governed by the coupling of 

the flow within the two-dimensional membrane to the hydrodynamics 
of the surrounding bulk fluid. We present a numerical scheme which 
is capable of solving this flow problem for arbitrarily shaped vesicles 
using the Oseen tensor formalism. For the particular problem of simple 
^ \ shear flow, stationary shapes are found for a large range of parameters. 

The dependence of the orientation of the vesicle and the membrane 
velocity on shear rate and vesicle volume can be understood from a 
simplified model. 



^ ! PACS numbers: 47.55.-Dz, 68.10.-m, 87.45.-k 

Vesicles are closed lipid membranes suspended in aqueous solution. If 
the composition of this solution inside and outside the vesicle is identical, 
it affects the equilibrium properties of the vesicle only by osmotically fixing 
the enclosed volume. Minimizing the bending energy of the membrane under 
the constraints of fixed enclosed volume and membrane area then yields the 
equilibrium vesicle shape at rest|I[]. The dynamics of this shape, however, is 
governed by the coupling of the flow within the two-dimensional incompress- 
ible fluid membrane to the hydrodynamics of the bulk fluid. Any theory of 
vesicle dynamics is complicated by the fact, that the boundary conditions for 
the three-dimensional Navier-Stokes equations have to be evaluated at the 
vesicle surface, which is moving with the fluid and whose shape is not known 
a priori. 



1 



Considerable practical interest in this problem arises from a more com- 
plex system, the red blood cell. The deformation of erythrocytes in hydro- 
dynamic flow fields is used as a measure of changes in the elastic properties 
of pathologically modified cells || The understanding of the dynamics 
of a single cell is also a prerequisite for the understanding of the rheology 
of blood. Since the erythrocyte has a complicated structure consisting of a 
fluid lipid/protein bilayer and the spectrin network, theoretical analysis has 
employed expansions for small flow rate [|j] or has focused on simpler model 
systems such as shells with stretching elasticity and liquid drops with homo- 
geneous surface tension ^, [7], ^. The dynamics of fluid membranes with 
fixed area and bending rigidity has only been studied in the quasiplanar |§, 



quasispherical [TU| and cylindrical [i f] geometry without external flow fields. 
Further analytical work includes the covariant treatment of the membrane 
incompressibility |T2] and the dynamic renormalization of material param- 
eters |U| in equilibrium. The problem for a closed vesicle with arbitrary 
shape and external flow is investigated in this Letter for the first time. 

For specifity, we investigate a linear shear flow f°(r) = 72 with shear rate 
7, as sketched in Fig. [I] (a), which can easily be realised experimentally @. 
The stationary state of a vesicle in such a flow resulting from our computation 
is shown in Fig. [1] (b). It is characterized by both a finite inclination angle 
9 between the longest axis of inertia of the vesicle and the flow direction, 
and a 'tank-treading' tangential motion of the membrane around the vesicle 




Figure 1: Vesicle in shear flow = 72. (a) Schematic drawing of the 
undisturbed flow, inclination angle 9, and coordinate and c. 

The coordinate axis pointing into the paper plane is y. (b) Stationary state 
obtained by numerical integration for v = 0.9, 7=1 s _1 , k = 10~ 19 J, 
77 = 1CT 3 Js / m -3 , and Rq = lO/xm. The velocities at the vertices are shown 
as arrows with length proportional to velocity. 
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with revolution frequency to. Both 9 and uj are found to depend strongly on 
the reduced volume v = V / (AnR^/?)) of the vesicle, where V is the enclosed 



volume and the surface area A determines the length-scale Ro = J A /An. The 
bending rigidity of the membrane k sets the energy scale |14] . Together with 
the bulk viscosity r), these quantities define an intrinsic timescale r = T]R^/ k. 
Thus, the dimensionless shear rate x = 7 r — IV^o/ K * s ^ ne secon d parameter 
of the problem. Figures |2| and |3| show the inclination angle 9 and the mean 
revolution frequency uJ [fl5[] , respectively, as a function of the reduced volume 



v for different values of x- 

We find numerically that for any non-zero shear rate small perturbations 
of the equilibrium shape do not allow stationary flow with the appropriate 
boundary conditions |16j . Thus, shear appears to be a singular perturbation 
that drives the vesicle towards a stationary non-equilibrium shape that can 
be distinctly different from the corresponding equilibrium shape. Further 
increasing the shear rate x only marginally affects shape and orientation 
angle 9. This behavior is remarkably different from the behavior of both 
erythrocytes, which go through a transition from tumbling to tanktreading 
behavior with increasing 7 H, and liquid drops, which deform to more and 



more elongated shapes as a function of 7 until they break up [|17j]. 

The stationary state was found to be independent of initial conditions. In 
particular, for v ^ 0.75, oblate discocytes are locally stable in the absence of 
an external flow field |Tj| but still transform to the same elongated shape as 



prolate vesicles when suspended in shear flow. We have tested this behavior 
for all v > 0.52, i.e. the entire range where the equilibrium discocytes do not 
self-intersect. 

Hydrodynamic calculations on typical length and time scales of mem- 
branes and vesicles are simplified by employing the Stokes approximation 
which is valid in the limit of small Reynolds number Re |[L9| . For a vesicle 
suspended in a shear flow with typical velocity u = 7-Ro> this number can be 
estimated by Re = puR Q /r] = p^R\jr]. Inserting the vesicle size R ~ 10/im 
and the viscosity rj = 10 _3 Js/m 3 and density p = 10 3 kg/m 3 of water, the 
limit of small Re corresponds to 7 10 4 s _1 . Thus, the hydrodynamic equa- 
tions for the velocity field v(r) of an incompressible Newtonian bulk fluid 
with external force density k(r) reduce to the instantaneous force balance 

r/V 2 v = Vp(r) + k(r). (1) 

The pressure field p(r) can be eliminated by using the incompressibility con- 
dition V • v(r) = 0. Due to the linearity of (P, v(r) can be written as a 
superposition of a background flow v°(r) and additional flows v l (r) arising 
from the interactions of the fluid with the membrane. 



3 




1 1 1 1 1 1 1 1 1 1 1 

0.5 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9 0.95 1 

V 



Figure 2: Inclination angle # between the x-axis and the longest axis of 
inertia as a function of reduced volume v for different reduced shear rates 
X =1 (A), 5 (0)> 10 50 (x), 100 (+). Numerical errors are smaller 
than the symbol sizes. The continuous line follows from Eq. (|8|). 
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Figure 3: Scaled revolution frequency ujt/x averaged over all vertices as a 
function of reduced volume v for different reduced shear rates x (Symbols 
are as in Fig. g). The continuous line is computed using the solution of (0). 
For volumes v < 0.7, a full revolution in the stationary state could not be 
computed due to constraints on computation time. 
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We describe the membrane surface by a triangulation with vertex po- 
sitions R a moving with the adjacent bulk fluid due to non-slip boundary 
conditions and impermeability of the membrane for bulk flow. Forces on the 
membrane are computed on these vertices from the energy 

g({R a }) = Y,MH 2 ) a + Y, a ° Aa ( 2 ) 

a a 

of the discretized membrane. The energy has two contributions. First, the 
discretized squared mean curvature (H 2 ) a ^Of describes the bending energy 



21] , [TjJ with bending rigidity k. Second, a locally varying isotropic tension 
a a which is conjugate to the area A a of the neighboring triangles around each 
vertex is needed to ensure local incompressibility of the membrane. 
The force K^ 3 at vertex R/ 3 then reads 



K 13 = 



<%?({R a }) 



and the force density is given by k(r) = J2p ~ R /3 )K /3 . 

The disturbance flow due to the presence of the vesicle can be computed 
using the Oseen-tensor O describing the flow field induced by a discrete point 
force K^ 3 in an infinite fluid medium. The total velocity at the vertex position 
R a then reads 

v(R Q ) = v°(R a ) + °( Ra > R-^K' 3 . (4) 
The matrix elements of the Oseen-tensor are p3[ 



with d = R Q — R^ and d = |d|. Here, the viscosity of the fluid rj is assumed 
to be equal for the fluids inside and outside the vesicle [EE . 

After dealing with the energetics of the membrane in the force bal- 
ance we have to treat the dynamics of the lipid membrane correctly. Two- 
dimensional incompressibility of the membrane leads to the condition that 
the area around every vertex remains constant under the dynamics, i.e. = 
d t A a = E /3 (<9A a /9R /3 )v(R /3 ). Using Eq. (|) for the velocities and Eq. (|) 
for the vertex forces, this condition yields an inhomogeneous system of linear 
equations for the local tensions a a . After the tensions are thus determined, 
moving the vertices by At v(r Q ) with a short timestep At respects the Stokes 
equations and all boundary conditions and thus gives a correct integration 
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scheme for the hydro dynamic problem which automatically conserves the 
volume of the enclosed fluid. Fluidity of the membrane is ensured by bond 
flips attempted in regular intervals [JZ7|]. 

This procedure turns out to be stable in the sense that (i) an initially 
smooth surface remains smooth during time evolution, and (ii) the distri- 
bution of triangle sizes and angles remains approximately constant in the 
stationary state. As a test, we relax surfaces at v° = and arrive at the 
known axisymmetric equilibrium shapes with an error of 0.2 percent for 
the energy and 1.5 percent for the (uniform) tension using a discretization 
with N = 337 vertices. Area and volume remain constant with an error 
below 0.1 percent for the longest runs. 

The numerical results for shape, orientation and membrane velocity in the 
stationary state may be understood using a simplified model. The numerical 
data show an approximate mirror symmetry of the induced tension o with 
respect to the plane defined by the axis denoted by a in Fig. [1| (a) and the 
y-axis. The axis perpendicular to this symmetry plane is denoted by c. We 
now assume that this symmetry holds exactly. The velocity field can be 
written as a superposition of the undisturbed shear flow v° and the flow 
v* induced by all forces on the membrane. The latter flow must have the 
mirror symmetry of the inducing forces. In a stationary state, the sum of 
these two contributions has to be tangential to the membrane surface. We 
evaluate this condition for the contour of the vesicle in the x-z-plane at y — 0. 
With s denoting the arclength of this contour, the contour is described by a 
vector (a(s), 0, c(s)) with the local tangential velocity @(s). We obtain four 
equations 

±/3{s)^j^- = vl + 7 cos 9(a(s) sin 9±c(s) cos 9), (6) 
ds 

/3(s)<^0 = — ^ sin 6> (a(s) sin (9 =b c(s) cos 6>) , (7) 

where the different signs apply to the upper and to the reflected point on the 
lower part of the vesicle, respectively. 

We can now eliminate the unknown quantity v l from (§) - (^) and obtain 
differential equations for da(s)/ds and dc(s)/ds. These equations describe 
an ellipse 

fl2 ( g ) C2 ( S ) 1 /ON 

r z cos z V r z sin 6 
If one now assumes axisymmetry of the entire shape, the length scale 
r and the angle 9 are uniquely determined by the constraints on area and 
volume of the vesicle. The shear rate 7 scales the velocity (3{s) and does not 
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influence the shape or orientation. We obtain the differential equation 



i28n 



d(3(s) 7 a(s) c(s) sin 2 9 cos 2 9 (cos 2 9 - sin 2 9) 
ds = a 2 {s) sin 4 fl + c 2 (s) cos 4 9 ( ' 

for the velocity along the contour by forming the derivative of the geometrical 
condition (d a(s) / ds) 2 + (d c(s) / ds) 2 = 1 with respect to s. This equation can 
be solved numerically, giving the revolution frequency uj = 2tt / § (3{s)^ 1 ds. 

In the spherical limit v ~ 1, Eq. (|8|) gives 9 « 7r/4 and r ^ Rq. Eq. (|9]) 
simplifies to /3(s) = 7-Ro/2 = const, which is equivalent to to = j/2 = x/2r. 
Thus, all limit values for a spherical vesicle are identical to the results for a 



rigid sphere or a fluid drop with infinite surface tension in shear flow [IS]. 

As Figs. H and § show, our simplified model yields good quantitative 
agreement with the results of the hydrodynamic integration. The results for 
different shear rate x collapse as expected. The spherical limit v ~ 1 can 
not be reached numerically, as the tension in the membrane diverges when 
exterior forces are applied to a sphere. For reduced volumes v iS 0.8, the 
symmetry assumptions in our simplified model are less justified. 

The results of the numerical integration contain information about the 
distribution of velocity and tension within the membrane not available from 
the simplified model. The revolution frequency UJ only measures the mean 
velocity. The local variations of membrane velocity show two remarkable 
features: (i) The velocity varies along the contour in a way that velocities 
are smaller towards the poles of the vesicle. This behavior can be understood 
qualitatively by formulating (|B|) and (|?p for a contour rotated out of the x-z 
plane, (ii) Fig. |] shows the lateral (y-)variation in the revolution frequencies 
u evaluated for individual vertices. Even though we have neither included an 
explicit shear rigidity nor shear viscosity, we find numerically effective shear 
in the membrane only for large deviations from an elliptical shape at small 
v. 

In the stationary state, the membrane tension a a is constant only for 
v° = 0. In a shear flow, we find that the largest values of the (negative) 
tension are in the middle of the vesicle, i.e. the vesicle responds to being 
pulled apart by the shear flow. For large shear rates bending rigidity 

is irrelevant, and the mean tension iV _1 J2 a °" a increases linearly with x, thus 
dominating the elastic energy of the membrane. 

In conclusion, we have developed a stable numerical scheme for calculating 
the time evolution of vesicles with arbitrary shape suspended in arbitrary 
flow fields. In the case of simple shear flow, we find an inclination angle 9 
independent of shear rate and a revolution frequency u linearly proportional 
to shear rate. Our approach breaks down for shear rates smaller than x ~ 
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Figure 4: Revolution frequency u (as a function of lateral coordinate y 
for reduced volumes v = 0.7, 0.8, 0.9, 0.975 for % = 10, r = 10s. The y- 
independence of uj for large v implies that adjacent membrane segments are 
again close together after a full revolution. 

k-QT/n ~ 0.04, where typical velocities of rotational diffusion /c B T/?7i?Q are 
comparable to the velocity of the shear flow. In this regime, one should 
incorporate thermal fluctuations. 

We thank T.M. Fischer for helpful discussions. M.K. and W.W. thank 
K. Kehr from IFF, Forschungszentrum Jiilich for hospitality. U.S. benefitted 
from interaction with S. Langer and M. Wortis on a related problem. 
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